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Abstract 

In the context of bilayer graphene we use the simple gauge model of Jackiw 
and Pi to construct its numerical solutions in powers of the bias poten- 
tial V according to a general scheme due to Kravchenko. Next, using 
this numerical solutions, we develop the Ermakov-Lewis approach for 
the same model. This leads us to numerical calculations of the Lewis- 
Riesenfeld phases that could be of forthcoming experimental interest for 
bilayer graphene. We also present a generalization of the Ioffe-Korsch 
nonlinear Darboux transformation. 



PACS numbers: 02.30.Hq, 11.30.Pb, 81.05.Uw 

1 Introduction 



The recent discovery by Novoselov et al [1] that a single atomic layer of graphite 
usually known as graphene can be stable under some special experimental con- 
ditions opened widely the realm of a soon-to-come technology of graphene-based 
electronic devices. The main transport mechanism of the charged carriers is in 
this case an enhanced quantum Hall effect accompanied by a non-zero Berry's ge- 
ometric phase that has been experimentally observed [2] . On the other hand, the 
quantum Hall effect and the Berry phase have been also discussed for graphene 
bilayers [3] that technologically could be even more interesting than graphene 
itself. Motivated by these works, in this paper, we introduce the more general 
Lewis-Riesenfeld phases for the case of a graphene bilayer. These phases are 
defined in a mathematical formalism called the Ermakov-Lewis approach, which 
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is shortly presented. For the bilaycr, we will make usage of the recent gauged 
Dirac model of Jackiw and Pi [4]. With the aid of a mathematical procedure 
due to Kravchenko [5] we present the results in a series expansion in the bias 
potential V applied to the bilayer graphene for any gauge potential A. We give 
detailed numerical results for a particular case of A. 



2 Jackiw-Pi model for bilayer graphene 

Recently Jackiw and Pi introduced a simple gauged Dirac model of bilayer 
graphene that can be reduced to the following coupled system of first-order 
differential equations [4] 

(A- - $(r))u(r) = Vv(r), (1) 
(D r + $(r))«(r) = -Vu(r), (2) 

where D r — , V is a constant representing the external bias voltage and 

*(r)=*/?-MW. (3) 
r r 

Physically function $ has two components: 

• a scalar field ip = that describes a particular dimerization called Kekule 
distortion [6] which can be also interpreted as characterizing the condensate 
arising from states bound by interlayer Coulomb forces between particles in one 
layer and holes in the other; the fingerprint of tp in <!> is the parameter k which 
should be an odd integer number. 

• a vector potential, A, which was first introduced by Jackiw and Pi for a 
monolayer graphene [7] to unpin the vortices and then studied in the context 
of bilayer graphene as well [4]. A is a gauge field satisfying some appropriate 
conditions in the origin and at infinity. In [4] it was supposed that 

A(0) = andyl(oo) = i. (4) 

It is easy to see that the spinor components u and v are necessarily solu- 
tions of the following second-order differential equations obtained directly from 
system (l)-(2) 

-(D r + <f>(r)){D r -$(r))u = V 2 u, (5) 
-(D r -<f>(r))(D r + <f>(r))v = V 2 v. (6) 

This couple of equations can be rewritten in the following way 

-D 2 r u + qi (r, V)u = 0, where qi{r, V) = $ 2 + - V 2 , (7) 
-D 2 v + q 2 (r, V)v = 0, where q 2 (r, V) = $ 2 - $' - V 2 . (8) 
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The latter form of the uncoupled system shows that the two functions u and v in 
the Jackiw-Pi model are supersymmetric partners in supersymmetric quantum 
mechanics and therefore many known supersymmetric results can be directly 
applied. 



3 Kravchenko's representation for the general 
Schrodinger solution 



Consider the Sturm-Liouville equation 



D r (pD r u) + qu = lo 2 u, 



(9) 



where p and q are complex- valued functions of a real variable r £ [a,b], p e 
C 1 (a, b) is bounded and non- vanishing on [a, b] and ui is an arbitrary complex 
constant. Suppose that there exists a solution g of the equation 

D r (pD r g ) +qg = 

on (a, b) such that go 6 C 2 (a, b) together with are bounded on [a, b]. 
Kravchenko [5] proved that the general solution of (9) has the form 



u = c\u\ + C2U2 
where c\ , C2 are arbitrary complex constants and 



ui = go V —rX 

^— ' 91. 



(n) 



U2 = 9o 



even n— odd n— 1 

with and X^ being defined by the following recursive relations 

1, 



n/;x(«- 1 )(0ff 2 (0« 



for an odd n 
for an even n 
for an odd n 

for an even n. 



(10) 



(11) 



(12) 



(13) 



nf a X^\09 2 (m 

This iterative scheme is appropriate for an easy implementation of numerical 
solutions and therefore for not exactly-solvable problems. 
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4 Application to the Jackiw-Pi uncoupled sys- 
tem 

To apply Kravchcnko's procedure to the Jackiw-Pi uncoupled system we start 
with the equation for u. It is easy to find a particular solution go from the 
factorized form of the u equation (5). We immediately get 

go' k (I A/ ^\ _ _h - k j ^ dr 



and the general solution u has the form u = C\U\ + C 2 u 2 , with 

— E V"' E V ■ (14) 

even n— odd n— 1 

where and are as in (12) and (13), respectively, for p = —1. 

The general solution v(r) of equation (6) can be obtained now by means of 
the Darboux transformation v = CiVi + C 2 u 2 = y9oD r gQ 1 (CiUi + C2U2) 

1 ~ jf(n) (fl) ~ X (n) 

By construction, «i and w 2 satisfy the following initial conditions 

tti(o) = 1, u[(a) = <&(a), 1*2(0) = 0, u' 2 (a) = 1, (16) 

with their Wronskian being equal to W — 1. Consequently, the functions «i and 
v 2 take the initial values 

«i(a)=0, «i(a) = -V, «2(a) = l ^(a) = (17) 
If one takes V = in (5), the formulas (14) give 



«i (r) = ci.g (r) and u 2 (r) = c 2 .g (r) / 

</ a 



(0 : 



where the expression in w 2 represents the well known construction of a second 
linearly independent solution of a Schrodingcr equation. 
When r — > 00, .4 — the system (l)-(2) reads as 

L> r u - Vu = 0, D r v + Vu = 0. 

A solution go m this case is an arbitrary constant, e.g., go = 1. Calculating 
and with the aid of (12)-(13) we obtain 
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{(—1) 2 r n , for an odd n ( (-1) ! r", for an odd n 

, X (n \r) = I 
(— l)f r™, for an even n \ (— l)%r n , for an even n. 

Now substitution into (14) and (15) gives us the following solutions (cf [4]) 
mi = cos Vr, u 2 = sin Vr, v\ — — sin Vr, v 2 = cos Vr. 

When r — > 0, and consequently ^4^0, the superpotential $ takes the form 
$ = Jp. For system (1), (2) with this particular potential, in [4] the following 
solution is presented 

u = VfJ(k-i)/2(Vr), v = -VrJ(k+i)/2(Vr), (18) 

which is easily obtained from (14) and (15) choosing go — r% and calculating 
m and v\. Indeed, for X^ we obtain the formulas 



X {n \r) = < 



("^lifeS^", for an odd n 

(-^^r", for an even n 



(fc+n-l)!! 

leading us to the solution (18) multiplied by a constant. 

We have now all the elements to proceed to specific calculations. Let us 
consider the Dirac system (l)-(2) choosing in (3) the vector potential A(r) = 
2 — \ e ~^ r which is in agreement with conditions (4). The parameter /j, -1 plays 
the role of a screening length. The superpotential (3) then becomes 4>(r) = 
^r-e~^ r and in this case the system (7)-(8) has the following form 

-Dlu + A e -*«Y* -/«■ _ I _ M ) w = V 2 (1Q) 

r 2r 2r r 

-D 2 r v + A e -fr(A e -f»- + - + n)v = V 2 v, (20) 
Ir Zr r 

The plots of u and v calculated within Kravchenko's scheme are presented in 
Fig.l. 



5 




5 The Ermakov-Lewis approach 

We will now expound on the Ermakov-Lewis (EL) formalism for the supersym- 
metric system (7)-(8). If the functions u and v are written in the form [8] 

u(r) = ap u (r) sm((p u (r, V) - 0), v{r) = jp v (r) sin((p v (r, V) - S), (21) 

where a, /3, 7and 5 are arbitrary constants, then both the amplitude functions p u 
and p v and the phase functions tp u and (p v are directly related to the Ermakov- 
Milne-Pinney nonlinear equations [9] 

Pu,v( r ) + <liA r > V)Pu,v{r) = 3 A ,, , (22) 

Pu,v V ) 

where A is an arbitrary constant, in the following way. While p u and p v are 
solutions of (22) with a respective potential, the phase functions ip u and <p v are 
expressed in terms of p u and p v as follows 

r 

V) = /^p) *'• CO) 

a 

For p u and p v there is an elegant representation [10] in terms of pairs of 
linearly independent solutions (1*1,1*2) and (i>i,t>2) of (7) and (8) respectively. 
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Namely, general solutions of (22) can be written in the following way 



p u = y Au\ + Bui + 2Cu lU 2 , p v = \J Avf + Bv\ + 2Cv Y v 2 , (24) 

where A, B and C are constants related by the equality AB — C 2 = with 
W being the Wronskian Uiu' 2 — u 2 u' 1 = V\v' 2 — v 2 v'y (which is a constant that 
does not change through the Darboux transformation). 

The solution of the Ermakov equations allows us to obtain another important 
dynamic invariant quantity I u and I v respectively, called the Lewis invariant 
[11]. I u is constructed in terms of any solution u of (7) and the corresponding 
solution p u of the Ermakov-Milne-Pinncy nonlinear equation (22) (and analo- 
gously I v is constructed in terms of v and p v ) 

Iu = \i^pJ + iPuU ' - ' Iv = l{^pJ + {PvV ' ■ (25) 



For the fixed interval [a, b] the expressions ^ip UfV (b,E) give us the quantum 
number function depending on the energy E [8], [12], [13] 



o 

N(E) = -[ -^-dr, 
ttJ (P{r) 



(26) 



used in standard quantum mechanics to identify bound states since at a bound 
state energy E n the function N{E) takes integer values N(E n ) = n + 1, n = 

0. 1.... and this does not depend on the choice of the constants A, B, and C. 
However, if the energy is different from E n then (26) is not uniquely determined, 

1. e., for E n < E < E n+ \ the function N(E) depends on the boundary conditions 
imposed to p. 

An important point of the EL formalism is that the phase function (23) also 
known as Lewis- Riesenfeld phase [14], has two components: the usual dynamic 
phase and a geometric phase that are given by the following formulas ([15]) 



A^(r) 



1 1 d 



pl iV (r') 2dr 



and 



( i V 
7 \Pu,v(r )- 



(r')\ / dp UjV (r') 



dr' 



dr' 



dr' 
(27) 



,sdp u ,v(r')\ ( dp UyV (r') 



dr' 



dr' 



dr', (28) 
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Any of these phases can be calculated by employing Kravchcnko's iterative series 
(14) and (15) for the p functions (24). Denoting 



E E 

even n—Oeven i—0 
oo n — 1 

E E 

odd n—leven i—0 

we get 



z!(n — z)! 

^•(i)^(n-i) 



-V 



i!(n — i)\ 



-V 



zu= E E 



odd n—lodd i—1 



i\(n — i)\ 



Pu(r) = g (r 

Introducing the notations 

X(i)X(n-i 

~* i\(n — i)\ 

i—1 



A ~ 2C ~ 

-J7T Zu + Bg 2 (a)Z u -—Z 
9o \ a ) v 



(29) 



= E E 

odd n—lodd i— 

oo n—1 

*«= E E 



■v n+ \ z v 



E E 

even n—Oeven 



JfWjf (n-i) 

i!(n-i)! 



v 



(i)jf(n-i) 

z!(n — z)! 



ocZtZ n— lev en i—0 

the amplitude function p v is calculated as follows 

, Bgl{a) 



Pv{r) = 



1 



9o(r) 



A 



V 



z 2 -^z 



(30) 



The particular values of the functions p u and p v can be obtained by relating the 
constants A, B and C to (16) and (17), that is 

p u (a) = VA, pja) = -j= , 



PvW 



— , PvW 



VB$(a) VC 

V 7b' 



(31) 



Based on (29) and (30), we have calculated the Lewis- Riesenfeld angular inte- 
grals (27) and (28) for the same decaying exponential profile of *4(r).The results 
for the geometric Lewis-Riesenfeld phases are plotted in Fig. 2 and Fig. 3 of 
the paper for A = 1 and A = 4, respectively, and for the same finite interval. 

One can check by direct calculation that the values of the Ermakov-Lewis 
invariant (25) do not change under tthe supersymmetric transformation, that is 
lu = Iv Indeed, noticing that I u and are r-independent quantities, we can 
calculate them at one point only, e.g., at r = a and with arbitrary constants 
d,C 2 and A, B, C. We have 

I u (a) = i (ACl + BC\ - 2CC1C2) , 

h{a) = K^r + BCl 2CCl ° 2 + 

Taking into account that C 2 — AB — A we obtain I u (a) = I v (a). 
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6 Generalization of the Ioffe-Korsch intertwin- 
ing formulas 

If the functions u and v are solutions of a coupled first-order system of the Dirac 
type then the corresponding amplitude functions p u and p v are related through 
nonlinear Dirac-Ermakov relationships that in a particular case A = 1 have been 
obtained previously by Ioffe and Korsch [13]. Consider the equality 

p 2 u = Au\ + Bu\ + 2Cuiu 2 . 

Substitution of u\ = —y(D r + 3>)ui and u 2 = —y(D r + <J>)v 2 leads to the 
relationship 

v 2 pi = $ 2 p 2 v + $(p 2 v y + s v , 

where S v = A (v^) 2 + B (v' 2 ) 2 + 2Cv[v' 2 . Taking into account (21), the functions 
v\ and v 2 can be expressed as v\ — ")ip v sin(<p„ — 5\) and v 2 = j 2 p v sin(ip v — S 2 ). 
Calculating the constants ji and Si according to (17) and (31) we get 

vi = -VBp v smtp v , 

v 2 = -^L (C sin ip v + cos ip v ) . 
v B 

Now taking the derivatives and substituting them into the expression for S v we 
obtain 

s v = (p' v f + — + (A - 1) I p' v sin ip v H cos ip v J , 

Pv \ Pv J 

where one can notice that (p' v sin ip v + j- cos ipSj = -g- (v[) 2 . 

Employing the same procedure with respect to p 2 u and taking into consider- 
ation that 

ui = (-C sin (p u + cos ip u ) , 
VA 

u 2 = VApu sin ip u , 

we arrive at the nonlinear Dirac-Ermakov system involving the amplitudes p u 
and p v and the phases <p u , ip v of the Dirac spinors u and v 

V 2 pI = {{D r + <5>)p v ) 2 + + (p' v sm<p v + — cos<p v ) , (32) 

Pv \ Pv J 



V 2 p 2 v = ((D r -<f>)p u y + — + (A-1) p; S in^„ + — cos^J , (33) 

Pu \ Pu J 

where again it is worth noticing that ( p' u si\np u + — cos ip u ) = -j (u 2 ) . Thus, 

there is a cross contribution to the squares of the intertwined amplitude func- 
tions of one-unit jumps in A which is given by the squares of the derivatives of 
the corresponding eigenfunctions. For A = 1 the system (32)-(33) reduces to 
the result in [13]. 



9 



7 Application to the Jackiw-Pi model for screened- 
exponential form of A(r) 

To associate the spinor components u and v to the amplitude functions p u and 
p v the relationship between the constants C\ and C2 and A, B and C should 
be fixed. For instance, choosing p u (a) — u(a), p' u {a) = u'(a) we get A = Cf, 
B = + Cf, C = C1C2. 

The corresponding A<^°* = (n + l/2)<p u<v and A^° m for A = 1 and A = 4 
are displayed in Fig. 2, Fig. 3. and Fig. 4. In addition, the amplitude functions 
for the same values of the parameters are shown in Fig. 5. Notice that since the 
problem is undefined at the origin because of the ~ singularity of $(r) the plots 
are displayed for the interval [0.l7r,7r]. For all plots we use the following three 
sets of values of the constants A, B and C: 

1. A = 1, B = A, C = 0. These values imply d = 1, C 2 = 0, 
u(r) — tii, i>(r) = Wi. 

2. A = ,4^, S = (A + l)$(a), C = -1. In this case Ci = -7=, 

C 2 = -$(a), ^(a) = 0. 

3. A = B = V 2 , C = leading to C x = ^g., C 2 = 0. 

The values of V, p and k are chosen as follows: V = 2.08236, p = 4j and fe = 1. 




Fig. 2. A^ ot (red) and A^* ot (blue) are 
plotted for three cases of the parameters A, B 

and C : — is for A = 1, B = A, C = 0, is 

for A = ^y, S = (A + l)$(a), C = -1, and -.- 

is for A = y$,B = V 2 , C = . Each of the 
cases is shown for A = 1 and A = 4. 



'V 

e' 

S 




Fig. 3. Plots of Aifl e om (red) and A^f om 
(blue) given by (28) for the same three cases of 
the parameters A, B and C as in Fig. 2. Each 
of the cases is shown for A = 1 . 
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Fig. 4. Plots of AifSf " 1 (red) and A^ eom Fig.5. The amplitude functions p u (red) and 
(blue) given by (28) for the same three cases of Pv (blue) given by formulas (29) and (30) for 
the parameters A, B and C as in Fig.2. Each the same parameters A, B and C* as in the 
of the cases is shown for A = 4. previous plots (A = 1). 



8 Summary 

After testing Kravchenko's numerical solutions in the context of Jackiw and 
Pi model of bilayer graphene we have developed an amplitude-phase approach 
of Ermakov-Lewis type for the same model and numerical solutions and for 
any radial profile of the gauge potential A. Next, we have chosen a Debye- 
screened form of the profile of A to illustrate our results. This choice could be 
of physical relevance when one considers the dynamic screening processes in the 
scattering of the charge carriers in the background of vortex configurations or 
other topological defects in the bilayer graphene. 

Our technique is general and therefore can be applied to any other physi- 
cal problem of similar mathematical structure. The geometric phase of Lewis- 
Riesenfeld type, which is more general than the Berry phase, is calculated for 
the first time in the graphene context. In addition, we provided a generaliza- 
tion of the Ioffe-Korsch nonlinear intertwining of the Milne-Pinney amplitude 
functions. 
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